Resonantly enhanced and diminished strong-field gravitational-wave fluxes 
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The inspiral of a stellar mass (1 — 100 M©) compact body into a massive (10 5 — 10 7 Mq) black 
hole has been a focus of much effort, both for the promise of such systems as astrophysical sources of 
gravitational waves, and because they are a clean limit of the general relativistic two-body problem. 
Our understanding of this problem has advanced significantly in recent years, with much progress in 
modeling the "self force" arising from the small body's interaction with its own spacetime deforma- 
tion. Recent work has shown that this self interaction is especially interesting when the frequencies 
associated with the orbit's 9 and r motions are in an integer ratio: Qg/Q r — /3g//3 r , with fie and 
fir both integers. In this paper, we show that key aspects of the self interaction for such "reso- 
nant" orbits can be understood with a relatively simple Teukolsky-equation-based calculation of 
gravitational-wave fluxes. We show that fluxes from resonant orbits depend on the relative phase 
of radial and angular motions. The purpose of this paper is to illustrate in simple terms how this 
phase dependence arises using tools that are good for strong-field orbits, and to present a first study 
of how strongly the fluxes vary as a function of this phase and other orbital parameters. Future 
work will use the full dissipative self force to examine resonant and near resonant strong-field effects 
in greater depth. 

PACS numbers: 04.30.-w, 04.25.Nx, 04.70.-s 



I. INTRODUCTION 



Our understanding of the two-body problem in gen- 
eral relativity has advanced substantially in the past 
decade. Besides the celebrated breakthroughs in numer- 
ical relativity 1-3] which have opened the field of bi- 
nary phenomenology in general relativity, there has been 
great progress in understanding the extreme mass ratio 
limit of this problem, when one member of the binary is 
much smaller than the other. This limit is of great in- 
terest in describing astrophysical extreme mass ratio bi- 
naries (a particularly interesting source for space-based 
gravitational- wave measurements) Q, and as a limiting 
form of the more generic two-body problem 0, Q . 

Most efforts to model extreme mass ratio binaries have 
focused on the computation of self forces (see Ref. Q 
for a recent comprehensive review). Consider a small 
body orbiting a black hole. At zeroth order in the small 
body's mass, its motion is described as a geodesic of the 
black hole spacetime. At first order in this mass, the 
black hole's spacetime is slightly deformed. This defor- 
mation changes the trajectory that the small body fol- 
lows, pushing it away from the background spacetime's 
geodesic. It is useful to regard the change to the tra- 
jectory as arising from a self force which modifies the 
geodesic equations typically used to describe black hole 
orbits. Conceptually, it is useful to split the self force 
into two pieces: A time-symmetric conservative piece, 
and a time-asymmetric dissipative piece. On average, 
the impact of the conservative contribution is to shift or- 
bital frequencies away from their geodesic values. The 



dissipative self force is equivalent, on average, to a slow 
evolution of the otherwise conserved constants (e.g., the 
orbital energy and angular momentum) which character- 
ize geodesic orbits. It makes the largest contribution to 
an orbit's phase evolution. The conservative piece makes 
a smaller (though still significant) contribution which ac- 
cumulates secularly over many orbits Q. 

Recent work by Flanagan and Hinderer Q (hereafter 
FH) using a post-Newtonian (pN) approximation to the 
self force together with fully relativistic orbital dynamics 
has shown that a small body's self interaction becomes 
particularly important near resonances. The background 
geodesic motion can be characterized by three orbital fre- 
quencies with respect to Boyer-Lindquist time: A radial 
frequency Q r , a polar frequency fig, and an axial fre- 
quency fi^. In the weak-field (large separation) limit, 
these three frequencies asymptote to the Newtonian Ke- 
pler frequency. In the strong-field, these frequencies can 
differ significantly, with Q r < fig < 0^. Resonant or- 
bits are ones for which the radial and angular motions 
become commensurate: flg/fl r = fig/fir, where fig and 
fir are integers with no common factors. On such or- 
bits, components of the self interaction which normally 
"average away" when examined over a full orbital period 
instead combine coherently, substantially changing their 
impact on the system's evolution. 

For the purpose of our background discussion, it is use- 
ful to include more details from FH's analysis of how reso- 
nant effects arise. Consider a body of mass \i moving on a 
bound trajectory near a Kerr black hole of mass M, with 
fx -C M, FH note that one can describe the motion of this 
body using action-angle variables and correctly account- 
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ing for how the integrals which parameterize geodesic 
orbits evolve due to the self force. Writing the angle 
variables q a = {qt,qr,<l9,q<p) (which describe motions in 
the t, r, 9, and (f> directions of Boyer-Lindquist coordi- 
nates) , and writing the integrals associated with geodesic 
motion Ji = (E,L Z ,Q) (with E the energy, L z the ax- 
ial angular momentum, and Q the Carter constant), the 
equations of motion describing the system are fioj 



have), or else because of other effects which resonances 



dg a 
dr 
dJj 

~d7 



L0 a (3) + egi 1 \q r ,qg,J)+O(e 2 ) , (1.1) 



eG\ l \q r ,qe,J) + 0(e 2 ) 



(1.2) 



The time parameter r is proper time along the orbit; the 
parameter e = /i/M, the system's mass ratio. The io r ,g.4> 
are fundamental frequencies with respect to proper time 
associated with bound Kerr geodesic orbits. The forc- 
ing functions g^P and G^ arise from the first-order self 
force. FH also include discussion of second-order forcing 
functions, which we do not need for this synopsis; see 
Ref. Q for further discussion. 

At order e°, Eqs. (|1.1|) and (|1.2|) simply describe 
geodesies of Kerr black holes: The integrals of the motion 
are constant, and each angle variable evolves according to 
its associated frequency. The leading adiabatic dissipa- 
tive correction to this motion can be found by dropping 



the forcing term g^ 1 and replacing G\ L> by (Gj iJ ), the av- 
erage of this forcing term over the 2-torus parameterized 
by qg and q r [13]. To compute this torus-averaged self 
force, it is sufficient to use the radiative approximation 
Q , which includes only the radiative contributions to the 
self interaction and neglects conservative contributions. 
The conservative contributions influence the motion only 
beyond the leading adiabatic order [H, [1(3, EH ■ 

Important post-adiabatic effects can be found by con- 
tinuing to neglect ga \ but now integrating Eq. (|1.2[) us- 
ing G] rather than its averaged variant. FH show that 



(i) 



for "most" orbits, G- lj is given by (G\ x> ) plus a rapidly 
oscillating contribution. Over the timescales associated 
with inspiral, this rapidly oscillating piece averages away 
and has little effect. The effect of the forcing term G^ 

is dominated by (G^ ) for "most" orbits. 

However, for some orbits, this averaging fails. When 
Qg/Q, r = /3g/f3 r , where fig and /3 r are small integers with 

no common factors, the contributions beyond (G^) are 
not rapidly oscillating. Such "resonances" instead can 
importantly modify how the integrals of motion evolve 
during an inspiral. A given binary is very likely to evolve 
through several such resonances en route to the final 
merger of the smaller body with the large black hole. 
A complete quantitative understanding of these resonant 
effects will thus be quite important for making accurate 
inspiral models. Prior to FH's analysis, several other pa- 
pers argued that such resonances may play an im por tant 
role in the radiative evolution of binary systems [l2t [l3j 
(albeit without quantifying the detailed impact they can 
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have on the evolution of a dynamical system [14 1 . 

Orbits in which Clg/Cl r take on a small-integer ratio 
have been studied in great detail by Grossman, Levin, 
and Perez-Giz [ll| , who called them "periodic" orbits and 
provided a fairly simple scheme for classifying their fea- 
tures. They have also demonstrated that periodic orbits 
may play an important role in laying out a nearly optimal 
computational strategy for sampling the parameter space 
of large-mass ratio orbits more generally [l6| . Following 
Ref. [9( , we will call them "resonant" orbits, reflecting 
the fact that our main interest is in understanding how 
their periodic structure impacts the self interaction. 

As a binary evolves through a resonance, its self inter- 
action and thus its evolution are modified compared to 
what we would expect if the resonance were not taken 
into account. The details of how the self interaction is 
modified depend on the relative phase of the radial and 
angular motions as the orbit passes through resonance. 
Because of this, resonances enhance the dependence of a 
binary's orbital evolution on initial conditions. Let the 
phase variable xo define the value of the orbit's 9 angle 
at the moment it reaches periapsis (see Sec. Ill Al for more 
details). On resonance, two orbits which have the same 
energy E, the same axial angular momentum L z , and 
the same Carter constant Q will evolve differently if they 
have different values of xo- 

FH estimate @ that the shift to the orbital phase in- 
duced by these resonances can be several tens to ~ 10 2 
radians (as compared to an analysis which neglects the 
resonances). That there is such a large shift, and that this 
shift may depend on initial conditions, is potentially wor- 
risome. Resonances could significantly complicate our 
ability to construct models for measuring the waves from 
extreme mass ratio inspirals. On the other hand, the 
detailed behavior of a system as it evolves through reso- 
nances may offer an opportunity to study an interesting 
aspect of strong-field gravity, providing a new handle for 
strong-gravity phenomenology. 

The "several tens to ~ 10 2 radians" estimate is based 
on applying pN self force estimates to strong-field orbits, 
a regime where pN approximations are generally inaccu- 
rate. It is thus of great interest to estimate the impact of 
orbital resonances using strong-field methods. One can 
in fact compute the dissipative piece of the self force at 
leading order in the system's mass ratio. Techniques for 
doing so with scalar fields were presented in Ref. fl7| : 
generalizing to the gravitational dissipative self force is 
not terribly difficult [HI, [l3[ • This computation uses so- 
lutions of the frequency-domain Tcukolsky equation [l8| , 
which has been used in recent years most commonly to 
study "flux balancing" [l9l - [2l| . the limit in which the 

driving force G\ is replaced by its torus average (G^). 

A full analysis based on a strong-field dissipative self 
force is underway, and will be presented later. Our pur- 
pose here is to present a snapshot of the impact of res- 
onances computed using strong-field techniques. In par- 
ticular, we demonstrate how the phase-dependent effects 
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discovered by FH appear when an orbit is exactly on res- 
onance, and we examine the numerical magnitude of the 
shift in the evolution of the integrals of motion due to 
this effect. In the language of Eq. (|1.2[) , we compare the 
on-resonance rate of change of the integrals of motion 
with the rate we would obtain if we torus averaged the 
forcing term 

We begin this paper by briefly reviewing the behavior 
of Kerr geodesic orbits in Sec. [TTJ Much of this material 
has been presented elsewhere, so we leave out most de- 
tails, pointing the reader to appropriate references. Our 
main focus is to describe how to find and characterize res- 
onant orbits. We then describe how to compute radiation 
from Kerr orbits in Sec. IIIII We first briefly review the 
Teukolsky-cquation-bascd formalism we use (Sccs. HTTAl 
MI C|) . and then describe how key details are modified by 
orbital resonances in Sees. MI Dl and MI El We describe 
two complimentary approaches to computing fluxes on 
resonance. Although formally equivalent (as we prove 
in Appendix lA"f . their implementation is quite different. 
Having both methods at hand proved useful to us in our 
numerical study. One aspect of the on-resonance com- 
putation (the evolution of Carter's constant Q) is suffi- 
ciently complicated that all details of this calculation are 
given in Appendix [B] 

Our results are given in Sec. IIV1 We begin by examin- 
ing how fluxes from individual modes (i.e., harmonics of 
the orbital frequencies) behave as a function of the offset 
phase of the radial and angular motions, xo- We show 
that the amplitude of a given mode, and hence the rates 
of change of conserved quantities associated with that 
mode, can vary significantly with \o- For example, the 
flux of energy from an orbit can vary by factors of order 
unity as xo varies from to 27r. The rate of change of the 
orbit's Carter constant can even change sign as xo varies. 
The total flux from a given orbit is given, however, by 
adding fluxes from many modes. When many modes are 
combined, much of the variation washes away; we find 
variations of a fraction of a percent in most quantities 
after summation. The amount of this residual variation 
seems to depend most strongly upon the geometry of 
the orbit's (r, 9) motion on resonance, in particular the 
topology of the trace in the (r, 9) plane. Orbits whose mo- 
tion in (r, 9) have a simple topology with few trajectory 
crossings in the plane (e.g., the fle/fl r = 3/2 resonance) 
tend to have relatively large variation in the integrals 
of motion; orbits whose motion has a more complicated 
topology with many trajectory crossings show much less 
variation (e.g., the Q,g/VL r = 4/3 resonance). We argue 
that this can be explained in terms of how the orbital 
motion tends (or fails) to average away variations in the 
source-term to the Teukolsky equation. 

Understanding how these fluxes behave exactly on res- 
onance is only the first step in building a complete strong- 
field understanding of how resonances impact inspirals. 
To go further, it will be necessary to examine how dis- 
sipation behaves as the system evolves up to and then 
through an orbital resonance. This analysis is just be- 



ginning; we briefly outline the approach we are pursuing 
in Sec. El 

Throughout this paper, we use "relativist's units," set- 
ting G = 1 = c. 



II. KERR GEODESICS AND ORBITAL 
RESONANCES 

A. Brief summary of general characteristics 

We begin by reviewing geodesic orbits of Kerr black 
holes, with a focus on aspects of this motion particularly 
relevant to our analysis. In most textbooks [for example, 
Ref. [3, Eqs. (33.32a)-(33.32d)], Kerr geodesies for a 
massive body are described using equations of motion in 
the Boyer-Lindquist coordinates i, r, 9, and <f>: 



, d9 

dr 



[E(r 



-A [r 2 + (L z - aEf + Q] 
R(r) , 

Q - cot 2 9L\ - a 2 cos 2 0(1 - E 2 ) 



= 0(9) , 

= esc 2 9L 

= $(r,0) 

= E 



(r 2 +a 2 ) 



aE 



2\2 



-1 - 



(2.1) 



(2.2) 

a 2 L z 
A 

(2.3) 



A 

+aL z (\ - 
T(r,6) . 



a 2 sin 2 i 



r 2 + a 2 



(2.4) 



In these equations, r is proper time along the geodesic, 
S = r 2 +a 2 cos 2 9, and A = r 2 — 2Mr+a 2 . The quantities 
E and L z are the orbital energy and axial angular mo- 
mentum, normalized to the mass fi of the orbiting body, 
and Q is the orbit's Carter constant, normalized to fi 2 . 
These three quantities are conserved on any geodesic. 

Along with the coordinate time t and proper time r, 
it is often very useful to work using a time parameter A, 
defined by dX = dr/Y,. The geodesic equations parame- 
terized in this way are 



dr \ 

5a .I = RW - 



d\ 



$(r, 



dt 
dX 



T(r,0). 



(2.5) 



By using A as our orbital parameter, the r and 9 coordi- 
nate motions are completely separated from one another. 
Proper time r couples r and 9 by the factor S; the cou- 
pling with coordinate time t is even more complicated. 
The parameter A is often called "Mino time," following 
Mino's use of it to untangle these coordinate motions 11 
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We have found it useful for many of our studies to 
introduce the following reparameterization of r and 9: 

pM 

r=— -, cosd = cos0 ro cos(x + Xo) ■ (2.6) 

1 + e cos y) 

These transformations replace the variables r and 9 with 
secularly accumulating angles ip and \. As ip and x 
evolve from to 2-7T, r and 9 move through their full 
ranges of motion. We define x = ip = at A = 0. 

Notice that we include an offset phase xo for the an- 
gular motion. Wc could also include an offset phase ipo 
for the radial motion, as well as initial conditions (po and 
to for the <p and t coordinates. We choose our time origin 
such that t = when A = 0, which means to = 0. Wc 
likewise choose (po = 0. Changing (p is equivalent to ro- 
tating around the black hole's spin axis, and can have no 
effect on the flux of energy and angular momentum from 
the system (although it introduces a phase offset to the 
system's gravitational waves). 

Finally, we choose ipo = 0, which amounts to setting 
A = at a moment that the orbit passes through pe- 
riapsis, r = r pori = pM/(l + e). The offset phase xo 
thus sets the value of 9 at periapsis. Previous work (e.g., 
[l9| ) has typically used xo = as well. The parameter 
set (^o, Yo, (pp, tn) is equivalent to the set (AJj, Aq, (po, to) 
used in Ref. [17[. Following this reference, xo = will 
label the "fiducial geodesic." Wc will use it as a reference 
geodesic for some of the calculations in Sec. IIIII 

In their original form, Eqs. (|2.1[) - (|2.4|) . Kerr orbits are 
parameterized (up to initial conditions) by the three con- 
served constants E, L z , and Q. The reparameterization 
(|2.6[) maps those constants to parameters that describes 
an orbit's coordinate geometry: semi-latus rectum p, ec- 
centricity e, and minimum angle 9 m . These quantities 
are likewise conserved along a geodesic. Schmidt [23[ 
provides closed-form expressions for converting between 
(E,L Z ,Q) and (p, e,6> m ). Either the set (E 7 L Z ,Q) or 
(p, e, 9 rn ), plus the relative phase xo, completely specifics 
a geodesic for our purposes here. 



B. Orbital frequencies and resonances 

Each orbit has a unique set of frequencies describing 
its motions with respect to r, 9, and (p. The frequencies 

Qr,e,4> = ^/T r fi,4, (2.7) 

are conjugate to the periods 1 expressed in coordinate 
time t] the frequencies 

Y r ,e,0 = 27r/ h. r ,e,4> (2.8) 



1 Describing the periods using Boycr-Lindquist time t is a bit com- 
plicated; T r g ^ really describe an averaged notion of the periods. 
See Refs. [23l , |24|| for more detailed discussion. 



are conjugate to these periods in Mino time A. These two 
frequencies are related by a factor T which describes the 
average increase in t per unit A: 

Qr,e,<f> = T^g^/r . (2-9) 

Details of how to compute T given (E, L z , Q) or (p, e, 9 m ) 
are given in Ref. J24|. One could also define frequencies 
conjugate to proper time r (see, e.g., Ref. [23| and dis- 
cussion in Sec. HJ, but the ft and T frequencies will be 
sufficient for our purposes. 

As background to Sec. IIII[ we compare a typical orbit, 
for which the ratio flg/fl r is some irrational number, to a 
resonant orbit, for which D,g/fl r = (3g/(3 r , where (3$ and 
j3 r are small integers with no common factors. Figure 
Q] shows the motion of three orbits, projected into the 
(r, 9) plane. In all cases, we have chosen p = 3.2758, 
e = 0.7, 9 m = 70°; the motion is thus bound to the range 
I.93M < r < 10.9M, 70° < 9 < 110°. 

In the right-hand panel, we have set the spin parame- 
ter a = 0.95A/. For these orbital parameters, this orbit 
has n$/£l r = 2.0311.... This is not a resonant orbit; 
notice that the roughly nine radial periods shown here 
do not close. The orbital trace in this case ergodically 
fills the (r, 9) plane. In the left-hand panel, we have set 
a = 0.9M, which yields Qg/Cl r = 3 — these orbits are 
in a 3:1 resonance. The two traces shown in this panel 
correspond to different choices of xo- The green trace has 
Xo = (so that 9 = 9 m = 70° at periapsis), and the red 
trace has xo = 7r /2 (so that 9 = 90° at periapsis). Both 
traces show roughly nine complete radial periods. By 
their periodic nature, their motion traces out Lissajous 
figures: No matter how long we follow these orbits, they 
trace out a 1-dimcnsional trajectory in the (r, 9) plane. 

Note that the geometry of the traces in the left-hand 
panel varies significantly as xo is varied. The topology 
of these traces remains fixed, however: In all cases the 
trace oscillates three times in the angular direction as 
it completes a single radial oscillation. As emphasized 
by Grossman et al. [Hj], the topology of resonant orbits 
is uniquely determined by their orbital parameters, by 
virtue of the integers (3g and (3 r that determine their pe- 
riodicity. We show some evidence in Sec. IIVI that the 
topology of resonant orbits directly affects the strength 
of their resonance. Simple orbits, which do not cross 
themselves often and do not "cover" much of the allowed 
(r, 9) plane, show large variations in their radiated fluxes 
as the phase xo as varied; more complicated orbits, which 
cross themselves many times and come "close to" much 
of the allowed (r, 9) plane, do not show such large varia- 
tions. 



III. GRAVITATIONAL RADIATION FROM 
KERR ORBITS 

Here we describe in detail how we compute radiation 
from strong-field orbits, with an emphasis on how reso- 
nances modify the "usual" behavior. We begin in Sec. 
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FIG. 1: Left: Lissajous figures describing motion in the (r, 6) plane on a 3:1 orbital resonance (a = 0.9M, p = 3.2758A/, e = 0.7, 
0m = 70°). The green trace has 9 = 9 m at periapsis; red has = tt/2 at periapsis. Approximately nine radial cycles are used 
to generate these traces. Right: Ergodic motion of a "normal" orbit. The orbit's geometry is identical to that in the left-hand 
panel, but we have changed the black hole's spin to a — 0.95M; this changes the ratio of frequencies to £le/tt r = 2.0311 . . .. 
Again, roughly nine radial cycles are shown here. Given enough time, the blue trace would pass arbitrarily close to all points 
in 70° < 6 < 110°, 2M<r< 12M. 



IIII A I by briefly reviewing the Teukolsky equation and its 
solutions. This material has been presented at length in 
several other papers, so we only give a summary. Our 
goal is to provide just enough detail to understand how 
the situation changes on resonance. Section IIII Bl de- 
scribes how to compute fluxes of energy E and angular 
momentum L z from Teukolsky equation solutions, high- 
lighting how this calculation must be modified for res- 
onant orbits. The analogous calculation for the Carter 
constant calculation is sufficiently complicated that we 
present its details in Appcndix[B] with a summary in Sec. 
HITCl Finally, Sees. UTTD] and UTTE] present two different 
ways to compute on-resonant fluxes. These methods are 
equivalent to one another, although their computational 
implementations are quite different. 



A. The frequency-domain Teukolsky equation and 
its solutions 

Our computation of the small body's self interaction 
uses the Teukolsky equation [l8j]. This equation governs 
the radiative components to a Kerr black hole's spacetime 
curvature, tpo and -04, which arise due to some perturbing 
source or field. In the relevant limits, identities make it 
possible to obtain all information about the field iJjq from 
-04, and vice versa, so we need only focus on one. The 
field ip4 is particularly convenient for describing radiation 



at infinity. 

Teukolsky showed [HI that, imposing the Fourier and 
multipolar decomposition 

/oo 
du^Rlm^SlmuiOy^-^ , (3.1) 
-°° lm 

where p = — l/(r — ia cos 6), a master partial differential 
equation governing ip^ separates. The function Si mu (0) 
is a spin- weighted spheroidal harmonic; Ref. [25| presents 
techniques for computing it to high accuracy. The radial 
function is governed by 

(I^P) ~ V ^ R l™ = -Wr,Xo) ■ (3.2) 

Equation (|3.2p is the Teukolsky equation (although that 
name is also used for the PDE that governs ^4 before 
separating variables). Setting the right-hand side of (|3.2[) 
to zero, we construct a pair of homogeneous solutions, 
Rimui (which is regular on the event horizon) and Rf^^ 
(which is regular at infinity). See Ref. [19| (hereafter 
DH06) for detailed discussion of how we construct these 
solutions, as well as for the potential V(r) appearing in 
Eq. (|3.2[) . From these solutions, it is straightforward to 
build a Green's function which can then be integrated 
over the source Timu to construct a particular solution. 

The source 7~i muJ is sufficiently complicated that we will 
not write it out explicitly; see DH06 for details. It is built 
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from projections of the stress-energy tensor for a small 
body orbiting the black hole, 



T a /3 



fXU a Up 



E sinOdt/dr 



6[r-r (t))6[6-e o (t,xo)}6[<t>-(l>o(t)] > 

(3.3) 

where \x is the mass of the small body, and u a are the 
components of its orbital 4- velocity. The subscript "o" on 
the coordinates in the delta functions stands for "orbit," 
labeling the orbit's coordinates (as opposed to a general 
field point, which we leave without a subscript). 

Note that Ti mu is a frequency-domain quantity. Be- 
cause it arises from Kerr orbital motion, it only has sup- 
port at frequencies uj m kn = mO,^ + kflg + nVt r , and is 
non-zero only for r min < r < r max , 9 m < 6 < it - 6 m 
[where r rmin = p/(l + e), and r max = p/(l - e); see 
Eq. (|2.6p ]. Once fully constructed, Timw has terms in 
<5[?' — r a (t)] and its first two radial derivatives; sec Sec. 
Ill of DH06. 

To understand fluxes from this system, our interest is 
in Rimui{f) m the limits r — > oo and r — t r + (the event 
horizon). These limits will allow us to deduce how the 
orbit evolves due to radiation to infinity, and due to ra- 
diation absorbed by the hole. As r -> oo, the homoge- 
neous solution Ry^ (r) approaches (modulo a power-law 
scaling) an outgoing plane wave. Likewise, as r — > ?'+, 
the solution Rp mu {r) limits to an ingoing plane wave. 
The particular solution we construct by integrating the 
Green's function over the source then takes the form 



Rlmui {r) 



Til 



i m (Xo)RZl m (r) 
Uxo)i?S m (r) 



00, 

r+, 



where 



ZLJxo) =C* dr 

Jr + 



mu) ( r j XO ) 



A(r') 



M2 



(3.4) 



(3.5) 



and where ★ can stand for oo or H . The symbol C* is 
shorthand for a collection of constants whose value is not 
needed here; see Sec. Ill of DH06 for further discussion. 

Next insert Timm into Eq. (|3.5[) and perform the r in- 
tegral. The result is a Fourier transform: 

/oo 
dte^-^It m Jr o (t),9 o (t, X0 )} 
-OO 

/OO 
dAe i( w r-mT^)A x 
-OO 

JLJr o (X),d o (X, X0 )} ■ (3.6) 



The function If mu introduced on the first line of Eq. Q3.6P 
is built from 7z mtlJ ; see Eqs. (3.30)-(3.33) in DH06 and as- 
sociated text for detailed discussion. On the second line, 
we have changed the integration variable from coordinate 
time t to Mino time A, and defined 

J*m»(ro,0 o ) = If m „(r ,e )T(r ,e )x 

e i[uAt(r o ,0 o )]-mA<t>(r o ,9 o )] _ >j\ 



[In any place that we indicate a dependence on (r D , 8 ), 
please note that this is shorthand for [r (X), 9 (X, Xo)]-] 
The function J* mLU (r 0l d ) is just If muj (r ,d ) rewcighted 
by T(r Q , o ) [the right-hand side of the geodesic equation 
(pO]) ]. and with the factor e i("At-mA4>) included. The 
functions At(r o ,0 o ) and A<fi(r o ,0 o ) are oscillatory con- 
tributions to the t and <fi pieces of the orbit: 

i Q (A) = rA + A*[r o (A)A(A, X0 )] , (3.8) 
O (A) = T A + A0[r o (A),0 o (A, X o)] . (3.9) 

Both At and Acb oscillate at harmonics of Tg and T r ; 
see Ref. [24[ for detailed discussion. 

The function Jf muj ( r o, S ) gathers all the pieces of the 
integrand for Zf mu] that can be described as harmonics 
of Tg and T r . As such, it is useful to decompose it into 
these harmonics: 

JLJr ,6 ) = J2 J ^ mk n(Xo)e- l{kTa+nT ^ , (3.10) 

kn 

where 

Julmkn(Xo) = / \ 2 / d\ dX" 

l Z7r J JO JO 

e l{krsXe+nT ^JLJr (X r ),0 (X e ,xo)}- 

(3.11) 

We have here taken advantage of the fact that Mino time 
completely decouples the r and motions from one an- 
other. We imagine that these two coordinates depend 
separately on two different Mino-time variables, A r and 
\ e . and integrate over a full period of each time. See Ref. 
[24| for detailed discussion of this trick. 

Next, combine Eqs. (j3U), (JXT]), ([3~TU]) . and (|3~TT|) to 
find 

27]- 

kn 

= Z lmkn(Xo)S(uJ - U) mkn ) . (3.12) 

kn 

On the last line, we have taken advantage of the fact that 
the delta functions mean that the RHS only has support 
at uj = uJmkm and we have defined 



2tt 

^uihnkn(Xo) = — JZlmkniXo) 



(3.13) 



^ dX s dY 

27rr J J 

(3.14) 



and 



Z lmkn(Xo) — Z Z mkn lmkn(X0) ■ 



(3.15) 
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Throughout this synopsis, we have explicitly shown the 
dependence on the relative phase Xo- To account for its 
influence on the amplitudes, let us first define 



Zlmkn = Zlmkn (XO — 0) 



(3.16) 



In other words, amplitudes with a check mark " are com- 
puted using the fiducial geodesic. As shown in Sec. 8.4 
of Ref. [T3|, the effect of xo is to introduce a phase: 



7* fvnl — pi£mfcn(X0) y* 



kn ' 



(3.17) 



where 



Ukn(xo) = kr e \ e o +mA$[r min ,0(-\ e o )} 

-uj mkn Ai[r min ,6(-\ 6 )} , (3.18) 

where A(f> is A<f> for the fiducial geodesic (and likewise for 
At), and where \® = Xq(xo) is the value of X e at which 
6 = 9 m . It is given explicitly by Eq. (3.75) of Ref. [nj]. 
On the fiducial geodesic, Aq = 0, and £, m kn — 0, as it 
should. 



B. Non- resonant fluxes of energy and angular 
momentum 



Once the coefficients Zf mkn (xo) are computed, it is 
straightforward to compute the fluxes of energy and an- 
gular momentum to infinity and down the black hole's 
event horizon. In this section, we present in some de- 
tail the calculation of fluxes to infinity, in order to high- 
light the step of the calculation which changes for radia- 
tion from resonant orbits. Resonance changes the down- 
horizon fluxes in exactly the same way. Since the down- 
horizon calculation is more complicated, we do not go 
through it in detail, but just summarize the final result. 

Using Eq. (|3.1[) and the definitions which follow, we 
find that as r — > oo, 



^ = \ E e^^ZtknSlrnUOy^- 
Irnkn 

= - ^ ip4,lmkn ■ (3.19) 



Imkr 



Here, Si m kn(&) is the spheroidal harmonic Si mul (6) for 
ui = uJmkn- As r — » oo, ip4 —> (l/2)(h + — ih x ). This 
means that 



ih-. 



2 v 1p4,lmkn 
Irnkn mkn 



(3.20) 



A very useful tool for understanding the energy car- 
ried by gravitational waves is the Isaacson stress-energy 
tensor [261 . 



Tirad 



16tt 



(d^h + dph + + d^h x dph x ) . (3.21) 



The angle brackets mean that the quantity must be aver- 
aged over several wavelengths; see Ref. [26| and references 
therein for detailed discussion of the averaging procedure 
necessary for Tj^ to be a proper tensor to the relevant 
order. On the first line, denotes a covariant derivative 
with respect to the background; on the second line, we 
take advantage of the fact that this limits to the partial 
derivative as r — > oo. 

The energy flux, our focus here, is given by 



dE c 



dt 



lim r 2 / T^ d n k dn 



= lim r 2 / T™ d dSl , 



(3.22) 



where n k is a radially outward pointing normal vector, 
and the index k is restricted to spatial directions. 
Combining Eqs. (|3~20f - (pT22|) . we find 



dE c 



y y 

dt I , J 47TCJ m/c „UWfc'r l ' 

\L7nkn I'm' k'n' I 

(3.23) 

-04 is the complex conjugate of ■04. The sum over I is 
taken from 2 to oo; the sum over m from —I to /; and 
the sums over k and n are both taken from — oo to oo. 
Similar rules apply to the sum over primed indices. The 
angle brackets on the left-hand side mean that this rate 
of change is to be understood as one which is averaged 
over appropriate orbital timcscalcs. 

Consider now averaging the right-hand side over sev- 
eral wavelengths. Assuming that each frequency cj m kn * s 
distinct (an assumption that is only true when we are 
not on a resonance), then this averaging forces m = m! , 
k = k' , n = n'. Using the fact that 



Re 



</'4 



4,imfen'04,i'm'fe'n' 



dn 



Slmkn(9)Sl> m kn{9)dQ — <5«< 



(3.24) 



we find 



dE* 
dt 



E = E ■ ( 3 - 25 ) 

Imkn mkn Irnkn 



A similar calculation focusing on T^ d gives us the flux 
of axial angular momentum: 



dLl 



At I ^— ' 47TW 3 i 

Irnkn man 



/ j ^z^lmkn 
Imkn 



(3.26) 



Notice that the phase £, m kn does not appear in Eqs. (|3.25[) 
and (|3~26)) . 

The calculation of fluxes down the horizon is more 
complicated. The Isaacson tensor is not defined in a 
black hole's strong field. Instead, one must use the fact 
that the curvature perturbation due to the orbiting body 
exerts a shear on the generators of the horizon, which 
increases the black hole's surface area. By the first law 
of black hole dynamics, this in turn changes its mass and 
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angular momentum; see Refs. |27J,l28( for detailed discus- 
sion. Assuming flux balance, we can then read out the 
down-horizon fluxes: 

\ _ V"^ \ Z lmkn\ 2 = \^ 

Imkn mkn Irakn 



dt 



(3.27) 



dL 



H 



dt 



E 

imkn. 



m \ Z l?nkn\ 2 _ fH 
Otlmkn _ / t n z,lmkr, 

mkn Imkn 



(3.28) 

We refer the reader to Eq. (3.60) of DH06 for the down- 
horizon factor ai m kn- 



D. Radiation from resonant orbits I: Merging of 
amplitudes on resonance 

On resonance, Vlg /fig = O r //3 r = O, and so kVtg + 
nQ r = Nil, where N = k(3g + n/3 r . An infinite number 
of pairs (fc, n) are consistent with a given TV. For a given 
value of m, all pairs (k, n) satisfying k(3g + n/3 r = N will 
have mode frequency u> m kn = ^>mN = rnCl,/, + Nfl. 

Revisiting Eq. (|3 . 19[) . this means that only three in- 
dices are needed to describe the radiation on resonance, 
rather than four: 



1> 



res 
4 



i 



£ N (xo)SimN(e)e^-^ «') , (3.33) 

ImN 



C. Non-resonant rate of change of the Carter 
constant 

The Carter constant Q is the third integral which al- 
lows one to separate the equations of motion for orbits 
of Kerr black holes. Unlike the energy and axial angu- 
lar momentum, there is no simple formula describing the 
"flux" of Carter constant carried by radiation. However, 
one can formulate how Q changes due to radiative back- 
reaction. Taking into account only the dissipativc piece 
of the self force and averaging over very long times, Sago 
et al. [2l| (hereafter S06) showed that 



dQ°° 
dt 

dQ H 
dt 



= E 



7 H I 

Imkn I 



mkn 



kT e ) 



Imkn 



j 3 

mkn 



= E 



a lmkn\Z^ kn \ X 



2ttuj 

(£mkn 

+ kT 9 ) 



(3.29) 



imkn 



where 

c, 



m(cot 2 0)L 2 



2 

a LU r , 



)E 



■ kit 



(3.30) 



(3.31) 



It is interesting that the rate of change of Q can be fac- 
tored into quantities that are encoded in the distant ra- 
diation {Z?u n and Z^ kn ) and quantities that are local 



to the orbital worldline (C r , 



^mkr, 



and Tg). Using 



Eqs. p.25|) and (|3.27p . these results can be written 



dQ* 
dt 



2 E E l™kn X 

+ kTg) /b} r , 



Imkn 



(3.32) 

where * is either oo or H. We go through the Sago et 
al. calculation of (dQ / dt) in some detail in Appendix [B] 
in order to understand how to modify this result on an 
orbital resonance. 

Note that the rates of change (dE*Jdt), (dL*/dt), 
and (dQ* /dt) are equivalent for non-resonant orbits to 
the three components of the torus-averaged forcing term 
(G^) introduced in the introduction, albeit using coor- 
dinate time t rather than proper time r to parameterize 
the rate of change. This equivalence breaks down for 
resonant orbits. 



where 



7* 



N 



(Xo) 



7" 

"In 



kn t 



(3.34) 



(k,n) 



and where (k, n)jy denotes all pairs (k,n) which satisfy 
kf3g + n/3 r = N. In Eq. (|3.33p . the sums over I and m are 
exactly as before, and N is summed from — oo to oo. 

Equations (|3.33[) and (|3.34l) tell us that, as we enter a 
resonance, modes of ip^ which were distinct combine with 
one another: "lines" in the gravitational-wave spectrum 
merge. Each mode's contribution to the combined ampli- 
tude (|3.34[) is weighted by its phase £ m fcn(xo)- Revisiting 
the calculation of the fluxes using Eq. (|3.33[) rather than 
([3T9j) . we find 



dE° 
~~dt 



ixo] 



\ z iLn(xo)\ 2 

imN 4 ™»^ 
E E lmN{Xa) i 



ImN 



dE 



H 



dt 



-(Xo) 



l^v(xo)! 2 



ImN ^N 
E E hnN(Xo) , 



ImN 



dt 



ixo] 



H^miv(Xo)| S 



47TW 3 



E 

ImN 

E ^z?(roiv(Xo) , 



ImN 



dL 



H 



dt 



-(Xo] 



V- ™\Z% N (xo)\' 



ImN 



N 



E L ^lmN(Xo) 



(3.35) 



(3.36) 



(3.37) 



(3.38) 



ImN 



(The factor ai m N appearing here is the same as aimkn 
introduced earlier, but with ui m kn replaced by oj m N-) 
Thanks to the dependence of Z* mN on the relative phase 
Xo, the on-resonance fluxes likewise depend on this phase. 
In Appendix [Bl we show how the calculation of dQ/dt 
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is changed due to an orbital resonance. The result is 
'dQ c 



dt 



-(Xo) 



E \ Z hnN(Xo)\ 2 r 
7, 9. L~mN 



ImN mN 



Re[^ jv ( Xo )^ Af (xo)] 



ImN 



2lTUJ mN 



(3.39) 



dQ H , A m mN \zz N (x»)\ 2 r 
xo) J — 2^ ^^"^ mAr 



f/f 



ImN 



(xo)3^U(xo)] 



(mJV 



AT 



(3.40) 



The factor £ m Ar is the same as £ m kn with Lo m kn replaced 
by uj m N- We have introduced the modified amplitude 



(fc,n) N 



Noti ce th at y* mN (xo) is similar to 2* mAf (xo) [compare 
Eq. (|3.34[) ]. but with each term in the sum weighted by 
k. Equations (|3.39[) and (|3.40[) are used in the following 
section to study how the Carter constant's evolution is 
affected by an orbital resonance. 



E. Radiation from resonant orbits II: The 
constrained source integral of a resonant orbit 



The method described in Sec. IIIIDI builds the on- 
resonance amplitudes Zf mN (xo) from the amplitudes 
Zt , which are normally computed with frequency- 
domain Tcukolsky equation solvers, such as that de- 
scribed in DH06. The only modification is the need to 
compute the phase £ m fcn(xo)- 

One can also compute the on-resonant amplitudes by 
modifying the integral for the amplitudes Z\ mkn . Doing 
so, we compute Z* mN (xo) directly, without reference to 
the amplitudes Zf mkn . We begin this calculation by car- 
rying over without modification the computation of Sec. 

irfrxi up to E q . 



As before, we decompose J* muj into Tg and T r harmonics. 
However, we now take into account how these frequencies 
are related on a resonance: 



T* - j* -i(kT e +nT r )\ 

" Imuj / j "ullmkn 

kit 

El* -i(kp e +np r )T\ 
u ulmkn c ' 

kit 

= ^ ] <^ulmN ( 



-iNT\ 



(3.42) 



On the second line, we've used the resonance relation 
Tg/fo = T r /(3 r = T. We then use N = k(3g + n/3 r , 
and change notation slightly to distinguish the source 
amplitude J*i mkn from its on- resonance variant S^imN- 



The result, Eq. ()3.42j) . depends on only one fundamen- 
tal frequency, T. As such, our integral for J^imN i s taken 
over only a single time variable A: 



Jui mN (xo) = ^ / WTrfA JLMX),0 o (X,Xo)]e mTX . 

(3.43 ) 

Finally, by combining Eqs. (JHSJ), ([3J]), (j3~42]) . and (|3T43|l . 
we define 



(3-41) Zf mN ( X o) 



2n 



T 



7* 

'- l UJ, 



UmN ImN 



(Xo) 



(3.44) 



2tt/T 



^Jf roWmJV [r(A)^(A, X o)]e 



iJVTA 



(3.45) 



In combining these equations, we find a proportionality 
to 5{uj — oj m Ar), which forces the RHS to have support 
only at uj = <jJ m N- Although it may not be obvious, Eqs. 
p.34[) and (|3.45p are equivalent. We show this analyt- 
ically in Appendix [AJ and will demonstrate it numeri- 
cally in the following section. A conceptually attractive 
feature of Eq. (|3.45[) is that the integrand is only eval- 
uated at the coordinates (r, 9) which the on-resonancc 
orbit passes through. Changing xo changes the points 
(r, 9) at which the integrand has support. This is how 
the dependence on xo enters Zf mN in this calculation. 



N 



However, Eq. (|3.45|) can only be used for orbits that 
are exactly on resonance. Indeed, in any other case, the 
3-index amplitude Zi mN is not meaningful since the on- 
resonancc condition k(3g+n[3 r = N is not met. A suitable 
generalization of Eq. (|3.34j) for slightly off-resonance or- 
bits can be used to understand the behavior of ipi as one 
approaches and moves through a resonance. As such, the 
sum of phase- weighted amplitudes, Eq. (|3.34[) . is likely to 
be more useful for understanding the resonant self inter- 
action in full inspiral studies. In any case, we have found 
having two techniques for computing 2* mN (xo) to be 
very useful. The codes which implement these two formu- 
lae are quite different, so it is reassuring that their results 
arc in agreement. As discussed at the end of Appendix IB") 
it appears that the modified amplitude 3^miv(Xo) can also 
be computed with a one dimensional integral by propa- 
gating the operator (d9/d\)dg under the integral in Eq. 
()3.45j) . We have not yet tested this, though it would be 
a worthwhile exercise to do so. 
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IV. RESULTS: HOW RESONANCES IMPACT 
RADIATION 

A. Variation of modes with xo> an d comparison of 
two computational techniques 

Wc now discuss examples illustrating how wave ampli- 
tudes and fluxes are affected by orbital resonances. All 
of our results are computed using a version of the code 
described in DH06, modified to handle resonances. 

Begin with Fig. [5J which illustrates how Zff nN and 
■^hnN behave as functions of xo- For this example, we 
have put a = 0.9M, p = 8.7744M, e = 0.7, 9 m = 20° (for 
which f2g/Q r = 3/2), and we have chosen I = 4, m = 3, 
N = 7. In all panels, the green curves show Z* mN com- 
puted using Eq. (|3.34l) ; the red dots show the same quan- 
tity computed using Eq. p.45|) . The two methods agree 
to numerical accuracy (roughly 6 digits 2 ). All examples 
that we have examined show that Eqs. (|3.34[) and (|3.45[) 
agree perfectly (as we would expect from the calculation 
presented in Appendix lA]). Having both methods at hand 
was quite useful for debugging the on-resonance version 
of our code. 

Besides showing the excellent agreement between our 
methods of computing Zf mN , Fig. [2] also illustrates how 
Z* mN varies with xo- For this example, we find that 
\Z^ nN \ varies by about 25% from minimum to maximum, 
and I-Z^tvI varies by about 40%. The associated energy 
fluxes, which are proportional to the amplitude's mod- 
ulus squared, varies by about 55% and a factor of two, 
respectively. 

Figures [3] and |4] give two examples of the on- resonance 
rate of change of orbital constants. We show Ef mN and 
QimN f° r two orbits about a black hole with a = 0.9. Fig- 
ure [3] shows the I = 2, m = 2, N = — 5 mode computed 
for an orbit with p = 3.2758M, e = 0.7, 6 m = 70°; in 
this case, Q,g/Q. r = 3. Figure |4] shows the I = 5, m = —2, 
N = 11 mode for an orbit with p = 4.5322M, e = 0.3, 
9 m = 45°, for which O e /Q r = 2. 

In both cases, the flux of energy to infinity varies by 
a factor of about 3.1. This agreement is a coincidence. 
The down-horizon flux shows more variety, varying by a 
factor of about 6.8 for the 3:1 resonance, and by a factor 
of nearly 10 3 for the 2:1 case. (This large variation is 
because the flux comes close to zero at xo — 4.7.) The 
variation in Q^° 2 -5 is especially interesting for the 3:1 
resonance: It is negative over nearly half the span of XO; 
but is positive elsewhere. This behavior is unique to the 
on-resonance form of Qi m jq , and arises from the fact that 
it contains a term proportional to Re [-Z^m/v^im/v] • Be- 
cause the amplitudes Zi m pf and yi m N can have different 
phases, the behavior of Q* m jy can be more complicated 



than the behavior of the energy or angular momentum 
fluxes. Those fluxes are both proportional to \Z* mN \ 2 , 
and hence are positive or negative definite. 

The horizontal dashed lines in these figures gives the 
rate of change that would be found if the resonance were 
neglected. In other words, it shows the rate of change 
one would find by simply combining in quadrature all 
of the 4-index amplitudes Z* mkn which contribute to the 
relevant 3-index amplitude Z* mN (xo)- Its value is the 
average with respect to Xq of the resonant flux: 

rrrk, no-ICS _ ±V I TP* J\B 

2t Jo 

= ^ P T Ef mN (Xo)^d X o. (4.1) 
27T J dxo 

Recall that the parameter Aq, introduced in Eq. (|3.18l) . 
sets the value of A at which 9 = 9 m . An explicit expres- 
sion for the Jacobian dX^/dxa is given in Eq. (3.76) of 
Ref. [TtJ ■ It is not difficult to show that this result must 
hold 3 : combining Eqs. (|3.34[) and (|3.35[) . we have 

-^imAf(Xo) = "j 2 \^Wnkn\ 2 

(4.2) 

The first sum in this expression is, as usual, taken over all 
pairs (k, n)jv, as defined earlier. The second sum is taken 
over the pair of pairs (fc, n)jv and (k',n')N, with k 7^ k' , 
77 7^ n'. The first sum is exactly E^j^°~ lc " . Using Eq. 
(|3.18p . we see that on resonance, 

Cmfen — £mk'n' = {k — k )T Aq . (4-3) 

Hence this term averages to zero, demonstrating the va- 
lidity of Eq. (|4.ip . Similar results hold for all of the other 
rates of change we compute in this paper. 

These examples show that the flux carried in each 
mode can vary significantly as a function of xo- This 
shows that in principle resonances can have a strong im- 
pact on gravitational-wave fluxes. Notice, though, that 
the detailed dependence of each mode on xo varies quite 
a bit from mode to mode. It would not be surprising 
if much of the variation cancels out after summing over 
many modes. Wc examine this in the next section, check- 
ing to see how much flux variation remains when many 
modes are added. 



2 It is not difficult to do the calculations more accurately than 
this [29l43l| , but 6 digits of accuracy is good enough for this first 
strong- field examination of this effect. 



At one point in our analysis, preliminary results indicated that 
averages did not respect Eq. 14. H . Gabriel Perez-Giz insisted 
to one of us (SAH) that this must be an error. Indeed, these 
preliminary results were wrong. 



11 




FIG. 2: Comparison of two methods to compute the on-resonance amplitudes Z* mN . All panels correspond to radiation from 
an orbit with parameters p = 8.7744A/, e = 0.7, 9 m — 20°, about a black hole with spin a = 0.9M. For this orbit, ile/Q r = 3/2. 
We have chosen I = 4, m — 3, N = 7. Left panels show -Zj| 7 , right panels show -Z437; top panels show the real part, bottom 
panels the imaginary part. Green curves show the amplitude computed by the on-resonance merging of amplitudes discussed 
in Sec. IIIIDI Sec. MI El red dots show the amplitude computed using the constrained source integral presented in Sec. MI El 
The two methods agree to numerical accuracy (roughly 6 digits in this case). 



B. Sum over many modes: Variation of total flux 

We now examine the variation in total flux on reso- 
nant orbits, computing the sums (|3.35|) and (|3.36[) . Those 
sums are taken over an infinite number of modes, which 
we cannot do in a numerical calculation. We instead 
truncate the sum over index I at Z max = 6; for orbits with 
e = 0.3, we truncate the sum over N at iV max = 50, and 
truncate at N max = 100 for e = 0.7: 

(max I Nmax 

£*(*o) -EE E HnNtXo) ■ (4-4) 

1=2 m=-l Af=-AT max 

We have not performed a careful convergence analy- 
sis, but have found that increasing l max and iV max only 
changes our numerical results by an unimportant frac- 
tion for the orbits we have examined so far. We do not 
claim our accuracy to be good enough for "production" 
purposes, but claim it is good enough to illustrate the 
physics that we present here. 

Figure [5] shows one example of how, after summing over 
many modes, E* varies as a function of xo- We put a = 
0.9M, and choose an orbit with p = 5.48622M, e = 0.7, 
and 9 m = 70°, for which Clg/Q r = 3/2. The fractional 
variation in E* we find is much smaller than the variation 
we saw in individual modes: the summed flux to infinity 
varies by about 0.2%, and the down-horizon flux varies by 
about 6.7%. The down- horizon flux is much smaller than 



the flux to infinity, so the variations arc dominated by the 
behavior E°°. The behaviors of L*(xo) and Q*(xo) are 
qualitatively similar to E*(xo), so we do not show plots 
for those quantities. 

Tables U - IIVI present the fractional variation in E* , 
L*, and Q* for several orbits about a black hole with 
spin a = 0.9M. Within each table, we fix e and 9 m . We 
look at large and small eccentricity (e = 0.7 and e = 0.3), 
and large and small orbital inclination 4 (9 m = 20° and 
9 m = 70°). We then vary p to study radiation emission 
from four different resonances, 3:1, 2:1, 3:2, and 4:3. The 
fractional variation in a quantity X is defined as 

AX = l^maxl - |*min| ^ 

(|X max | + |X min |)/2 ' V ; 

where X max / min is the maximum or minimum value X 
takes as xo varies from to 2ir. 

The pcak-to-trough variation (|4.5I) in the fluxes is an 
important quantity that determines several properties of 
the resonances. First, the "kicks" in E, L z , and Q that 
occur as a system spirals through a resonance are di- 
rectly proportional to the variation (|4.5p [32J. Second, 
there are two qualitatively different types of resonances 



4 Note that smaller 9 m implies a more highly inclined orbit; 9 m = 
90° is an equatorial orbit. 
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a = 0.9M, p = 3.2758M, e = 0.7, B m = 70' a = 0.9M, p = 3.3758M, e = 0.7, 9 m = 70" 




Xq Xq 

FIG. 3: On-resonance variation of the rates of change of orbital energy (left panels) and Carter constant (right panels) in the 
I = 2, m = 2, N = -5 mode for an orbit with p = 3.2758M, e = 0.7, 6 m = 70°, a = 0.9M (for which Q. e /VL r = 3). Top 
panels give the flux to infinity, bottom ones give flux down the horizon. The dashed line in all panels shows the value that 
would be obtained if the resonance were neglected (i.e., simply adding in quadrature the various 4-index amplitudes Z* mkn that 
contribute to the 3-index amplitude Z* mN ). In all cases, the flux varies considerably with the phase \o- The variation in Q°° 
is especially interesting in this case, changing sign at xo — 0.28 and xo — 3.01. We do not show £*(xo) for this mode, since it 
is identical to E*(xo) modulo a factor of m/w m jv. 



that can occur in systems of this kind: a simple linear 
resonance in which the kicks depend sinusoidally on the 
phase parameter xo (cf. the final equation of FH), and a 
nonlinear variant in which the dynamics is rather more 
complicated. For the nonlinear scenario, it is possible 
to have a "sustained resonance" in which the system be- 
comes trapped near the resonance for an extended pe- 
riod of time [33l |34| . Our numerical results show that 
AX < 1 at least over all of the parameter space we have 
surveyed so far, which indicates that the resonances are 
always of the simple, linear kind. This agrees with post- 
Newtonian analyses [HJ . 

Some interesting trends are apparent from these tables. 
First, notice that in all cases the down-horizon variation 
is quite a bit larger than than the variation in the quan- 
tities to infinity. However, in all cases, the magnitude 
of the down-horizon fluxes is substantially smaller than 
the magnitude to infinity. The total variations are thus 
dominated by the fluxes to infinity, consistent with the 
results shown in Fig. [5j 

Second, notice that the largest variations are seen in 
either the 2:1 or 3:2 resonances (always the 3:2 resonance 
for orbits with e = 0.3, but either 3:2 or 2:1 depending 
on which quantity we examine for the orbits with e = 
0.7). The variations are consistently smallest for the 4:3 
resonance. This behavior correlates with the shape that 
a resonant orbit traces in the (r, 9) plane. Figure [5] shows 



these orbital tracks for the four orbits presented in Table 
[H For simplicity, we only show tracks for xo = ir/2. 

The contrasting shapes of the 2:1 and 3:2 orbits on one 
hand, and of the 4:3 orbit on the other, are particularly 
noteworthy. The 4:3 resonant orbit (bottom right) traces 
a rather complicated Lissajous figure which comes "close 
to" many of the (r, 9) points which are accessible given 
(p,e,9 m ). This complicated trajectory samples much of 
the accessible domain in r and 9. Appealing to the con- 
strained integral method of computing Z*i m ^ (cf. Sec. 
IIIIE[) . we can say that the motion effectively averages 
out the variations in the integrand by passing close to so 
many accessible points. 

By contrast, the trajectory for the 2:1 and 3:2 res- 
onances (top right and bottom left) are much simpler. 
These trajectories do not come as close to so many points 
in their allowed domain, and so do not average the varia- 
tions in their integrands as effectively The trajectory for 
the 3:1 (top left) resonance is similar to that for the 2:1 
case, but with an additional angular oscillation at small 
radius. This extra oscillation enhances the averaging as 
the orbit moves through a particularly strong-field part 
of its domain. Not too surprisingly, the flux variation in 
this case is generally intermediate to the others. 

Beyond the fact that orbits with simple shapes in the 
(r, 6) plane tend to show stronger resonances than or- 
bits with more complicated shapes, we do not as yet see 



13 



a = 0.9M, p = 4.5322M, e = 0.3, B m = 45° a = 0.9M, p = 4.5322M, e = 0.3, fl m = 45" 




FIG. 4: On-resonance variation of the rates of change of orbital energy (left panels) and Carter constant (right panels) in the 
I = 5, m = -2, N = 11 mode for an orbit with p = 4.5322M, e = 0.3, d m = 45°, a = 0.9M (for which Q.g/Q. r = 2). Top panels 
give the flux to infinity, bottom ones give flux down the horizon. The dashed line gives the value found when the resonance is 
neglected. As in Fig. [3] we see that Ei mN and Qi mN vary quite a bit as xo sweeps from to 2n, with minima near zero in this 
case for the down-horizon quantities. 
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TABLE I: Variation in flux for orbits with e = 0.7 and 9 m = 20° about a black hole with spin a = 0.9M. We vary p to examine 
a sequence of orbital resonances from Qg/Q r = 3 to £lg/Q r = 4/3. Columns 3-5 show the fractional variation in energy flux, 
axial angular momentum flux, and Carter constant rate of change arising from the down-hole fields; the fractional variation is 
defined precisely in the text. Columns 6-8 repeat this information for these fields at infinity, and columns 9-11 give the 
fractional variation for the totals (infinity plus horizon). The variations are largest for the 3:2 resonance and 2:1 resonances 
(depending on which quantity we examine), and smallest for the 4:3 resonance. 



strong evidence of any trend which would allow us to pre- 
dict which resonances will tend to be "strong" (i.e., ex- 
hibit large variation in orbital parameter evolution) and 
which "weak." Consider for example the rate of change 
of orbital energy, Ai? tot . As we go from high inclination 
to shallow and from high eccentricity to low, we see that 
AE tot goes from large to small: It takes the value 1.03% 
for high eccentricity, high inclination (Table IJ); 0.167% 
and 0.303% for the mixed cases (Tables ITT1 and UTT)) : and 
the value 0.131% for the case of small eccentricity, shal- 
low inclination (Tabic UV)) . This appears to suggest, at 
least roughly, that the strength of the resonance is cor- 
related with the degree of radial and angular motion. 

However, no such pattern is seen when we exam- 
ine AL* ot and AQ tot . For L z , the high eccentricity, 
high inclination case again produces the largest variation 



(0.489%, in the 3:2 resonance of Table H]). However, the 
low eccentricity, low inclination case produces the second 
largest variation (0.123%, in the 3:2 resonance of Table 
|1V[) . These values of e and 8 m likewise produce the largest 
and second-largest variations in the Carter constant (al- 
beit in different resonances). 

We do not yet have a compelling way to explain these 
trends (or lack of trends) in the resonances' strength, so 
we leave this mystery to future work. 



V. CONCLUDING DISCUSSION AND FUTURE 
WORK 

In this analysis, using a Teukolsky-equation-based for- 
malism good for exploring radiation produced by strong- 
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TABLE II: Variation in flux for orbits with e = 0.7 and 8 m — 70° about a black hole with spin a = 0.9M. As when e = 0.7 
and 9 m — 70°, the variations are largest for the 3:2 resonance and 2:1 resonances (depending on which quantity we examine), 
and smallest for the 4:3 resonance. 
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TABLE III: Variation in flux for orbits with e = 0.3 and 8 m = 20° about a black hole with spin a = 0.9M. In this case, the 
3:2 resonance shows larger variations than all other cases; the 2:1 resonance is surprisingly weak, given its strength in other 
examples we have seen. As usual, however, the 4:3 resonance shows the least amount of variation among all the resonances 
that we consider. 



a = 0.9M, p = 5.48622M, e = 0.7, 5 m = 70" 3:1 resonance 2:1 resonance 




FIG. 5: Variation of total energy flux, both to infinity 
(top) and down the horizon (bottom) for an orbit with 
p = 5.48622M, e = 0.7, 8 m = 70°, a = 0.9M (for which 
fig /fl r = 3/2). After summing over many modes, the varia- 
tion is significantly reduced: the flux to infinity only varies by 
about 0.127%, and that down the horizon varies by roughly 
1.6%. The variations in L1.(xo) and Q*(xo) are qualitatively 
similar, so we do not show them. See TablellTlfor more details. 



field orbits, we have confirmed the picture that on reso- 
nance the gravitational- wave driven evolution of a binary 
can depend strongly on the relative phase of radial and 
angular motions. A binary in which this relative phase 



FIG. 6: Trajectories in the (r, 8) plane for the orbits discussed 
in Tabled] We put xo = n/2 for these plots. The 4:3 resonance 
shows the smallest flux variation of those considered here, and 
has the most complicated trajectory. This orbit comes "close 
to" enough points in the (r, 8) plane that it averages over 
much of its accessible domain. By contrast, the 3:2 and 2:1 
orbits have simple trajectories and do not effectively average 
over their domain. Fluxes from these orbits tend to show 
the largest variation with xo- The 3:1 orbit is similar to the 
2:1 orbit, but with an additional angular oscillation at small 
radius which enhances orbital averaging. This orbit generally 
shows intermediate flux variation compared with the other 
cases. 
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TABLE IV: Variation in flux for orbits with e = 0.3 and 9 m = 70° about a black hole with spin a — 0.9M. The case is 
qualitatively similar to most of the others, with the 3:2 and 2:1 showing the largest degree of variation (depending on the 
quantity being examined), and the 4:3 case showing the least. 



has the value 7r/2 as the system enters resonance may 
evolve quite differently from an otherwise identical sys- 
tem in which this phase is 3tt/2 entering resonance. A 
typical extreme mass ratio binary can be expected to pass 
through several orbital resonances en route to its final co- 
alescence. That their evolution through each resonance 
depends strongly on an "accidental" phase parameter has 
the potential to complicate schemes for measuring grav- 
itational waves from these binaries. 

We find that the degree of variation depends strongly 
upon the topology of the orbital trajectory in the (r, 9) 
plane. Of the cases we have studied in detail, the or- 
bital plane trajectory of resonances like Clg/Q, r = 3/2 
have a simple topology. This trajectory does not cross 
itself very often, and does not come close to many points 
in the plane. Such resonances do not effectively average 
out the behavior of the source to the wave equation. As 
such, if the source varies significantly over an orbit, there 
can be a strong residue of this variation in the associated 
radiation. By contrast, the trajectory of resonances like 
Qg/Q r = 4/3 has a more complicated topology, crossing 
itself many times, and more completely "covering" the 
plane. In these cases, the orbit comes "close to" many of 
the allowed points in the (r, 9) plane, which quite effec- 
tively averages out the source's behavior. 

Although instructive and a nice validation of our abil- 
ity to examine resonances, these results are not enough 
to truly assess the importance that resonances have in 
a strong field analysis. We must be able to analyze a 
system as it evolves through a resonance, and thereby 
integrate the full "kicks" in the integrals of motion E, 
L z , and Q imparted to the system as it passes through 
resonance. To do this, we have begun expanding our 
Teukolsky code to compute, in the frequency domain, the 
instantaneous components of the dissipative or radiative 
piece of the self- force. Our formulation is based on the 
discussion in Refs. [infill, and is similar to the prescrip- 
tion described in Ref. [17| f° r scalar self forces 5 . This 
will allow us to study how a real inspiral is affected as 



One might be concerned about gauge ambiguities associated with 
the gravitational self force. As shown by Mino these ambi- 
guities disappear when one averages the self force's effects over an 
infinite time. In a two-timescale expansion [id ], such ambiguities 
remain, but are suppressed by the ratio of the timescales. 



we evolve through each resonance using results that are 
good deep in the strong field. The results shown in this 
paper are a first step toward this, demonstrating that our 
strong-field toolkit can be used to study resonant effects. 
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Appendix A: Proof: Equivalence of methods for 
computing on-resonant amplitudes 

In this appendix, we prove that Eq. p.45[) . the 1-D in- 
tegral for the on- resonance 3- index amplitude 2* mN (xo), 
is equivalent to Eq. p.34[) . the on-resonance amplitude 
expressed as a sum of 4-index amplitudes Zi mkn (xo), each 
of which is computed using the 2-D integral (|3.14|) . 

We begin with Eq. (|3.6[) . which we repeat here: 

POO 

J —oo 

(Al) 

Recall that the "o" subscript on r and 9 means that those 
are quantities along the orbit, and as such vary at har- 
monics of the frequencies T r and Tg. We can thus ex- 
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pand J* mul in a Fourier series: 

JU, = Y. J ^ kn (Xo)e-« kT ° +nr ^, (A2) 

kn 

= Y.^i mN (xo)e- lNTX . (A3) 



N 



Equation (|A2[) holds for arbitrary orbits. Equation (|A3[) 
only holds on resonance, when Tg = /3gT, T r = /3 r T. 

Because Eq. (|A2j) remains valid for resonant orbits, in 
the resonant case 



^Cimiv(Xo)e ?JVTA — ^ ^Zrnfcn(Xo) 



-i(fcT()+nT r U 



iV 



(A4) 

(The notation "=" means that this equation is true only 
on resonance.) Multiply both sides by e tN TA and inte- 
grate from to 27r/Y. On the left-hand side, we have 



Jo N 



2tt 

Y 



NN' 



N 



2lT ( \ 

~y~<JwlmN'\X0) 



(A5) 



To do this operation on the right-hand side, first note 
that by the resonance condition we must have 

kTg + nT r = (k/3 e + n/3 r )T . (A6) 

Using this, the integral for the right-hand side becomes 

f27r/T 



el [Af'-(fc/3 e +«/3 r )]TA rfA 



27T 



kn 



2jt 



(A7) 



j 



The notation (k, n) jy means that the sum is over all pairs 
(fc, n) which satisfy k(3g + n(3 r = N'. 



Next, use Eqs. ([335]) . (|5T3|) and (j3jH]), invoke Eq. 
(|3.34j) . drop the primes on the index N, and equate (|A5[) 
and ([AT]) . The result is 



(A8) 



(fc,n)iv 



which proves that the 1-D integral and the sum of 2-D 
integrals are equivalent for resonant orbits. 



Appendix B: Evolution of the Carter constant 



The third conserved quantity associated with orbits of 
Kerr black holes is the Carter constant, Q. Rearranging 
Eq. (j2~2j) . we write 



( ^ 



(Bl) 



Reference [21( (S06) first demonstrated how to compute 
the long-time-averaged evolution of Q, at least for non- 
resonant orbits. In this appendix, we revisit their calcu- 
lation in some detail in order to see clearly how it will 
have to be modified for resonant orbits (modifying some 
details to be in accord with our notation). We then ex- 
amine how the calculation changes when we are on an 
orbital resonance. 



1. Setup 



We begin with the first line of Eq. (3.18) of S06. It relates the averaged rate of change of the Carter constant, per 
unit Mino time, to the Kerr metric's Killing tensor K a/3 and to a radiative field ^rad which is constructed from the 
perturbation to the Kerr spacetime metric: 

/ — — \ = lim — — / d\— — = Van — — [ dX 

\ dX I l^oo 2L J_ L dX l^oo 2L J_ L 

We refer the reader to S06 for a detailed derivation of Eq. (|B2[) . and defer discussion of the radiative field Vl / la d(a;) to 
Sees. IB 3l and IB 4l The coordinate x represents a general spacetime field point; x — > z(X) means to take this general 
point to the orbit's worldline z(X). 

The other quantities appearing in Eq. (|B2[) arc as follows: First, K is a variant of the Carter constant, given by 



rati 



(B2) 



K = Q+ (L z - aE) 2 



(B3) 
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It is related to the Kerr metric 's Killing tensor by 

K = K aP u a u p , (B4) 

where 

K ali = 2SV Q m' 3) - a 2 cos 2 9g af) . (B5) 
The tensor g a @ is the Kerr metric, and m a are components of the Newman-Penrose tetrad leg, 



in 



iasine r e 1 , icsc6 

m = , m = —= : , nv = —= ; . (Bo) 



\/2(r + ia cos 9) ' V2{r + ia cos 9) \/2{r + iacosl 

Ovcrbar denotes complex conjugate. The quantity u a is the 4-velocity promoted to a spacetime field: 



(fit, u r , u e , u*) - (-B,±#)/A, ± y/etf), L z ) , (B7) 

where R(r) is defined in Eq. (|2.1j) . and <d(6) in Eq. (|2.2[) . Notice that our fie differs from that used in S06. This is 
due to a difference in the definition of the potential O (it describes motion in 9 here, but motion in cos# in S06). The 
field fi Q reduces exactly to the 4-velocity u a when we take the limit of the field point x to the worldline z(X). 



2. General simplification 



We now take the first steps in simplifying Eq. (|B2j) . These steps are the same for both resonant and non-resonant 
cases; we specialize to those cases in Sees. IB 3l and IB 41 
We begin by focusing on the integrand of Eq. (|B2|) : 



2Y l K ap u a d, 













x—>z{\) 



rail 



2£a 2 cos 2 9u a d l 



rad 



Use the fact that u a = u a in the limit x — > z(X), and that Y.u a = dx a /d\. Expanding trS a ra^\ we find 



2 (A) 



2Y>K a0 u a d B ( ) = 2£ 



[L z - asm 2 9E) ( esc 98^ + adt 



d9 
dX 



do 



rad 



2a 2 cos 2 6 



d 



dX 



vrad 



(B8) 



(B9) 



[For brevity, we omit x —> z(X) in Eqs. (|B9|) and (|B10[) . though it should be understood that this limit is taken.] The 
right-hand side of Eq. (|B9[) can be simplified significantly by combining the term in d9/dX with the final term: 



2A 
dX 



- 2a 2 cos 2 1 



d /* rad \ d9 n * rad dO n 

= 2—dey iad -2——deZ 
dX 2j dX 



dX 



-2a' 



dX 



cos 2 9- 



rad 



s 



2a 2 



rad 



(BIO) 



The third term on the right-hand side of Eq. (|B10|) is a total derivative in d/dX. Thanks to the periodic nature of all 
the relevant terms, it will not contribute to an averaging integral of the form (|B2[) . and may be discarded. Using 



= a 2 da cos 2 



dX 



a d0 x 2 

9 = —dg COS I 

dX 



(Bll) 



we see that the second and fourth terms on the right-hand side of (|B10I) cancel; only the term in c^e^rad remains. The 
integrand simplifies to 



2EK af3 u a d l 



1> 



rad 



[L x - a sin 2 9E) (esc 2 9d^ + ad, 



The radiative field ^rad can be broken into an "out" and a "down" component: 



^Tx de 



*rad 



!(A) 



(B12) 



(B13) 
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These two fields are in turn computed from mode functions $/ mw (discussed in more detail momentarily) as follows: 



*ss(x) = J ^Ei^kw / 



c.c. , 



(B14) 



down 
rad 



/ duj 4iui 2 Prl 



(B15) 



In Eq. (|B15|) . p m = cj — mQjj, where = a/2Mr+ is the angular velocity associated with the event horizon. The 
abbreviation "c.c." means complex conjugate. See S06 for further discussion and derivation of these forms of the 
fields ^°ad an d ^rad™- We largely focus on the "out" field, which is related to radiation at I + . Extension to the 
"down" field, related to radiation on the event horizon, is straightforward. 

To proceed, we use two equivalent forms for §^(3;) evaluated in the limit x — > z(X); both are described in more 
detail in S06. The first is up to a constant factor the complex conjugate of the integrand in the expression (|3.6[) for 
Z H ■ 



Irnuj 



(B16) 



Here T is the factor introduced in Sec. lIIBl that converts the mean accumulation of Mino time to the mean accumulation 
of coordinate time. Equation (|B16|) is Eq. (3.11) of S06, translated into our notation 6 ; the scalar-case version of this 
equation is Eq. (9.20) of Ref. [17T| . Using the Fourier series expansion (|3.10j) of J^ nu) , integrating with respect to A, 
and combining with the definitions (|3.13[) and (|3.15p gives 



[*(*')] 



n k 



(B17) 



and so 



*S3(x) = jdu 



Imkn 



x) 



c.c. 



A similar simplification describes ^rad™^)- Combining this with Eqs. (|B2|) and (|B12|) . we obtain 



(B18) 



dK°< 
~dX 



. ^ 2iuj 3 , 

\lmkn mkn 



d0 

(esc 2 9L Z - aE)dd, + a(L z - aE sin 2 0)d t + —d 9 

aX 



Lmkr 



(B19) 



where 



(IjOUt 



(The superscript "00" is because we focus on the "out" field.) 
We next manipulate the term in dg in Eq. (|B19[) . by invoking the second form for ^°^ u (x), which is 



(B20) 



The value of /zmfcn (?"j #) is not important for our purposes; see S06 [Eq. (3.20) and nearby text] for further details. 
We have changed notation from S06 slightly to highlight the fact that this function depends on Z, m, k, and n; this 
is important for generalizing to resonant orbits. We now evaluate on the worldline x —> z(A), and use the following 
explicit representations of the motions in t and <f>: 



t(A) = rA + At r (A) + Afy(A), 
0(A) = T A + A0 r (A) + A^(A), 



(B21) 



cf. Eqs. (|3.8p and (|3.9[) above and Sec. 3 of Ref. [l7j]. Here the function At r is periodic with period A r and Atg is 
periodic with period Ag, etc. This gives 



QZlnW = /;mfen[r(A), 0(X)} cxp {-zA(fcT e + nT r ) - iu mkn [At r (X) + At e (X)} + im[Acj> r {X) + A0 e (A)]} . (B22) 



Note that there are two errors in Eq. (3.11) of S06: the sign of 
the exponential is nipped, and the coefficients Z are of the usual 
type II3.15J I rather than the required more general type Il3.14j l 
with ui 7^ w m fe n . See Eq. (IB30t below. 



19 



We next define a mode function of two variables (A r , X 9 ) by 

<CL(A r , X e ) = fi, nkn [r(X r ), 9(X 9 )} exp {-ikT e X 9 - inT r X r - iuj mkn [At r {X r ) + At e (X 9 )} + im[A<p r (X r ) + A0 e (A 9 )]} . 

(B23) 

This function is determined uniquely by the following three properties: First, it reduces to the expression (|B22|) when 
evaluated at A r = X 9 — A; second, it is biperiodic, with a period of A r in A r , and of A e in X e ; and third, it is a 
continuous function of the geodesic's parameters. The first two properties are sufficient to guarantee uniqueness for 
non-resonant orbits, but not for resonant orbits since the different periodicities become degenerate. Adding the third 
property is sufficient to restore uniqueness for all orbits, since resonant orbits form a set of measure zero in the phase 
space. See Refs. [13, Hit f° r more details on the mapping between functions of A and functions of (A r , X e ). 

Next, differentiating the explicit expression (|B23|) with respect to X e , we obtain the following identity relating the 
differential operator d/dX and the partial derivative operators dg, dt and acting on $° m \ n '- 



d0„ d dAt 6f . dAfo 



(B24) 



We now use the identity (|B24[) to substitute for the {d6/dX)dg term in Eq. (|B19I) . This yields 

'dK c 



dX 



\ ^ Imkn 
\lmkn mkn 

esc 2 9L Z - aE 



dA<f) 6 
~d~X e ~ 



aL z — a E sin 



dAt 6 
~dX°~ 



d t + ikTe + -tt-t 



d 
dJ 9 



link i 



+ C.C 



(B25) 



Using Eqs. (3.43) and (3.58) of Ref. [17| it is not difficult to show that 

dAcj), 



esc 2 9L Z 



aE 



aL z - a 2 Esm 2 I 



dX 9 
dAt 



(esc 2 8L Z 



aE) 



esc ff)L, - aE 



dX 6 



= {aL z 



a 2 E sin 2 6) = aL z 



a z E (sm 1 



(B26) 
(B27) 



Combining this with Eq. (|B25P and using the replacements 3$ — > im, dt — > —ioJmkn gives 

dK c 



dX 



^Imkn 
\lmkn mkn 



iM m kn + ikT e + -fTn 



d_ 

dX 6 



imkn 



C.C. 



(B28) 



where we have defined 



Mmkn = m((csc 2 0)L Z - aE) - au mkn (L z - a£(sin 2 0)). 



(B29) 



Next, from Eqs. ([B16| . (j3~T0| . (|3~T3| and (f3~T5l) we obtain an expression for ^° m \ n {X). Extending this to a function 
of A r , A as above gives 



^ 1 _ 



JAkT e X" JAnT r X T 



Jmk+Ak.n+An^ 



(B30) 



An.Afc 



Combining this with Eq. (|B28[) yields the final result 



dK<* 
~dX 



7 H 

^ j^'Y* .1 Imkn r?H 



Imkn Afc, An 



LJ 3 
rnkn 



7 tl iAkTgX iAnT r \ T , 

^ui^nlmk+Akm+Ant e < 



(B31) 



Here it is understood that the averaging procedure is to first evaluate at A r = X 9 = X and then average over A. In 
Sec. IB 3[ we evaluate this average for non-resonant orbits, and reproduce the results of S06. In Sec. IB 41 we do so for 
a resonant orbit and find an appropriately modified variant of their formula. 
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3. Non-resonant result 

We evaluate the expression (|B31|) at A r = X = A and then evaluate the average over A defined by Eq. (|B2|) . The 
term labeled by An, Afc is proportional to 

lim -L f dAe iAfcT 6 A e iAnT r A = Um Si[( AfcT e + AnT r )L] , (B32) 

where Si(x) = sin(a;)/a;. Since the frequencies Yg and T r are incommensurate for non-resonant orbits, the combination 
AfcTg + AnT r will be nonvanishing for (Afc, An) ^ (0,0), and the right hand side will vanish. Thus the only non- 
vanishing term will be the term with An = Afc = 0. Another way to think about this is that we are averaging over a 
curve which is ergodically filling up the torus parameterized by A r and X s , and so the curve average can be replaced 
by an average over the torus, 

lim — / . . . dX -> -f-^r / / ... dX r dX e . (B33) 



l^oo 2L J_ L ' ' ' (2?r) 2 

Applying this torus average to the expression (|B31|) again forces An = Afc = 0. Now using the definition (|3.15|) wc 
obtain the final result 

« ) = T E '3 [Mmkn + kTg] + C.C. (B34) 
' Imkn nUJ mkn 

Because all the terms on the right-hand side of (|B34[) are real, the complex conjugate simplifies to an overall factor 
of two. We take the long-time average, so 

dK\ I dK\ 

*r)= r (*7- (B35) 

Further, by Eq. (|B3)l, 

dK dQ s / dE dL z \ , 

— = -± + 2{aE- L z ) a- -± . B36 

dt dt K ' \ dt dt J y ' 

Combining Eqs. (pT25)) , (|3~2"6]) . ([3"3Tj) . ([B29| together with Eqs. ([B"34]) , (|B35|) . and ([B36| . we finally obtain 

= E (m(cot 2 0)L a - a 2 w mfe „(cos 2 6)E + fcT 9 ) 

1 Imkn rnkn 

fpoo 

= 2 E — £«ta ■ (B37) 

UJ?nk7i 

Imkn 

A similar calculation focusing on the "down" modes yields 

(^r) = E aim * nlZ f knl2 {m(cot 2 e)L z -a 2 u lmkn (co S 2 e)E + kTe) 

O \ ^ Imkn r {J^QQ\ 

— ^ l-mkn ■ (OdO) 

, , ^mkn 
Imkn 

The factor ai m k n is introduced in Sec. lIIIB] on the second line, we have used Eqs. (|3.27[) and (|3.31[) . Equations (|B37[) 
and (|B38|) are the same (modulo minor changes in notation) as Eq. (3.26) of S06. 

4. Resonant Q 

We now return to the general formula (|B31[) evaluated at A r = X = X and compute the average over A for the case 
of resonant orbits. Before evaluating this average we first simplify the sums over Afc and An by rewriting them in 
terms of fc' = fc + Afc, n' = n + An. We also make the replacements 

E-E E > E-E E . CB39) 

kn N (fc,n)jv k'n' N' (k',n') N i 
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where the indicated sums are taken over k, n satisfying kf3g + n(3 r = N and over k', n' satisfying k'fy + n'f3 r = N'. We 
note that the quantities A4 m kn an d u mkn depend on k and n only through N, and write these as M. m N and <^> m N- 
Finally using the definition (|3.34j) of the amplitudes Z* mN , the expression (|B31|) reduces to 

^) = (^EE E [^ + fc%]^^ Wmfe ,„,e^V^ + cx\ (B40) 

Next we note that the argument of the exponential is 

i\{AkT g + AnT r ) = i\T(Ak/3 e + An/3 r ) = i\T(N' - N). (B41) 
Evaluating the average over A enforces N = N', and the result is 

^t) = £E E [M mN + k'T e ]^Z^ knlmk , n ,+c.c.. (B42) 

' ImN (k',n') N mN 

Now since to m kn — ^mN — ^mN'i the factor of Z^ imk'n' can ^ e simplified to Z^n^'n'- The expression (|B42|) can 
then simplified further by defining the new amplitude 

3Cv(xoH E kZ l H mkn ( X o)= J2 ke l ^^Z^ nkn . (B43) 

(fc,n) N (k,n) N 

Compare this with Eq. (|3.34p : yf^^ixo) is similar to Zi m ^(xo)^ but with each Z^ kn weighted by k. In terms of this 
new amplitude the result simplifies to 

/ dK°° \ r 

( ~a\ ) = E ^jr- [M mN \Z* nN {x Q )\ 2 + T z£ iW ( Xo )3Cv(xo)] + c.c. (B44) 

ImN m ^ 

Applying Eqs. f)B35|) and (|B36|) . we at last find the rate of change of Q for a resonant orbit: 

'^P) = E «^{A^&(xo)| 2 + T e Re[Z^(xo)&(xo)]} . ( B45 ) 
where 

C mN = m(cot 2 0)L Z - a 2 uj mN (cos 2 6)E. (B46) 
Repeating this exercise for the "down" modes yields 



^jr) = E {£ m Ar|Z°°( Xo )| 2 + T e Re [Z^( Xo )^( Xo )] } . (B47) 



ImN m 



It is interesting to compare our final result for the on- 
resonance evolution of Q, Eqs. (|B45|) and (|B47|) . with the 
equivalent results for the non-resonant case, Eqs. (|B37|) 
and (|B38|) . The first two terms in both expressions for 
(dQ/dt) are essentially the same; going from the non- 
resonant case to the resonant case is simply a matter of 
promoting the 4-index non-resonant amplitude Z^ mkn to 
the 3- index resonant amplitude Z^ mN . 

The final term in the two cases is quite different, how- 
ever. In the non-resonant case, the final term is propor- 
tional to kTg. In the resonant case, the index k can- 
not appear in the final result, which can only depend on 
the indices I, m, and N. This is accounted for in the 



definition of the amplitude y* mN , Eq. (|B43[) . In both 
the non-resonant and the resonant cases, this final term 
arises from the action of the operator (d6/dX)dg on the 
radiative field # ra d [sec Eq. (|B12[> ]. 

As Appendix [X] made clear, the 3-index amplitude 
Z£ mN can be computed directly as a 1-D integral, Eq. 
(|3.45[) . or can be computed as a sum of 4-index integrals, 
Eq. (|3.34[) . each of which is computed from the 2-D in- 
tegral (|3"3I]l . Our definition (|B43|) of yf mN is clearly 
analogous to Eq. (|3.34[) . writing this 3-indcx amplitude 
as a sum over 4-index amplitudes. 

Might it be possible to compute the 3-index amplitude 
directly, in a manner analogous to Eq. (|3.45[1 ? We be- 
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lieve the answer is yes: We simply need to propagate the We have not yet tested this, 
operator (dO/d\)dg under the integral sign in Eq. (|3.45p . 
In other words, we speculate that 

v T / >27r / T dQ 
yLAxo) = f J d\-d e JLJr(\),e(\ 7 xo)V NTX . 

(B48) 
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